#include<iostream>
#include<cmath>
#include<cstdlib>
#include<iomanip>
#include<fstream>
#include<sstream>
#include<string>

using namespace std; 
      
int main (int argc , char * argv [])
{



int total_number=8;
double g[15]={1e-6,1e-6,1e-6,1e-6,1e-6,1e-6,2e-6,3e-6,4e-6,5e-6,6e-6,7e-6,8e-6,9e-6,1e-5};
double diffu[15]={0.05,0.02,0.01,0.005,0.002,0.001,0.001,0.001,0.001,0.001,0.001,0.001,0.001,0.001,0.001};

double var1[11];
        
 	char poreFileName[128]="INPUT_CG_";
	char poreFileName2[128]="run_";
 ostringstream name;
for (int i=1;i<=total_number;i++)
{
        
        name.str("");
	name<<poreFileName<<i<<".inputdat";
	ofstream out;
	out.open(name.str().c_str());
	



//=========================================================================
//============================MC_CG========================================
//=========================================================================




out<<"========3D LATTICE BOLTZMANN MPI CODE--MULTI COMPONENT CG======="<<endl;
out<<"val_cap2.dat         :Geometry filenameEmily_Berea_y_614_764_7.dat"<<endl;
out<<"phase2_"<<i<<".dat           	:Initial components distribution"<<endl;
out<<"20 51 3               	:nx ny nz"<<endl;
out<<"1000000		     		:Maximum time step"<<endl;
out<<"7.0	             		:dx (um) Resolution (for Permeability calculation)"<<endl;
out<<"0	        	      	:Pressure Or Velocity Boundary (1=YES, 0=No)"<<endl;
out<<"0				:Psi constant BC"<<endl;
out<<"1.0e-6 0.0 0.0e-6	          	:body force for x,y,z"<<endl;
out<<"0 1.0 0 1.0027	    		:Pressure Boundary in X direction ()p=c_s^2*/rho"<<endl;
out<<"0 1.0 0 1.0	    		:Pressure Boundary in Y direction ()"<<endl;
out<<"0 1.0 0 1.0	    	:Pressure Boundary in Z direction (Format detials can be found within this file)"<<endl;
out<<"0 0.0 0 0.0	   	:Velocity Boundary in X direction (Format detials can be found within this file)"<<endl;
out<<"0 0.0 0 0.07	   	:Velocity Boundary in Y direction (Format detials can be found within this file)"<<endl;
out<<"0 0.0 0 0.07	   	:Velocity Boundary in Z direction (Format detials can be found within this file)"<<endl;
out<<"0 0			:Psi constant BC in X 0=OFF, 1=Ini psi, 2=Neibourghing,3=Rand"<<endl;
out<<"0 0				:Psi constant BC in Y"<<endl;
out<<"0 0				:Psi constant BC in Z"<<endl;
out<<"0.01    		       		:Viscosity (Component A, psi=1)"<<endl;
out<<"0.1    		       		:Viscosity (Component B, psi=-1)"<<endl;
out<<"-1.0     		      	:Contact Angle Cos(Theta) (Positive=1 wetting, Negative -1 wetting)"<<endl;
out<<"1.0e-2    		      	:Surface tension (Kappa)"<<endl;
out<<"0.0 0.0 0.0			:initial velocity for x,y,z"<<endl;
out<<"7.80177e+06	        	       	:Permeability (Single Phase mD)"<<endl;
out<<"=========OUTPUT==CONTROL==================="<<endl;
out<<"1		         	:Permeability writing (1=yes, 0=no)"<<endl;
out<<"1                       	:Direction of Permeability Calculation (1=X, 2=Y, 3=Z)"<<endl;
out<<"50		         	:Frequency of results writing (interval in time steps)"<<endl;
out<<"0               		:Memory Saving mode for output(0=No, 1=Yes,for large outputs, x,z opposite)"<<endl;
out<<"50000				:Freqency of velocity field writing (in vtk format, -1=no velocity writing)"<<endl;
out<<"50000				:Freqency of density field writing (in vtk format, -1=no density writing)"<<endl;
out<<"50000				:Freqency of concentration writing (in vtk format, -1=no density writing)"<<endl;
out<<"==============ADVANCE==PARAMETER=========="<<endl;
out<<"0 1.0 4.0	              :Self define lattice velocity: 0=DEFAULT,dx, dt ((u_x+u_y+u_z)*dt/dx<=1 Courant Number)"<<endl;
out<<"./"<<i<<"_	:OUTPUT PATH,DEFAULT: ./ (INCLUDE / AT THE END)"<<endl;
out<<"2            :PRESSURE AND VELOCITY BOUNDARY CONDITION OPTIONS: 0,1,2,3: EBC_S,EBC_D,TOLKE_BC,NEBC_D"<<endl;
out<<"0  500	:MULTI-COMPONENT STABALIZER: (a,b) a=0=OFF, a=1=ON, BODY FORCE APPLIED AFTER b steps"<<endl;
out<<"-1               :PRESET SATUATION, 0--1 initial component distri not needed, -1=OFF(for Component A,1)"<<endl;
out<<"0				:PRESET VALUE FOR BUFFET AREA, 0=NO,1=COMP A,-1=COMP B (valid when preset satuation)"<<endl;
out<<"0.0				:Relative permeability calcualtion 0..1 (psi>=value, cal the flux for Comp1)"<<endl;
out<<"0 0 0				:Pemeability calculation Partially  (1=ON, 0=OFF)"<<endl;
out<<"308 9				:Permeability calculation partially Starting point and Ending point in X"<<endl;
out<<"0 100				:Permeability calculation partially Starting point and Ending point in Y"<<endl;
out<<"0 100				:Permeability calculation partially Starting point and Ending point in Z"<<endl;
out<<"======================BACKUP==CONTROL===================="<<endl;
out<<"-1		                   :BACKUP FREQUENCY (-1=NO BACKUP,0=Backup at the end of computation)"<<endl;
out<<"0                              :INITIALIZATION WITH BACKUP DATA (0=OFF, 1=ON)"<<endl;
out<<"======================GEOMETRY READING======================="<<endl;
out<<"0		:Geometry Reading format, 0=decimal,1=binary"<<endl;
out<<"======================PRESSURE OR BODAY FORCE SETTING (CHANGE WITH TIME)============="<<endl;
out<<"0 1			:0=OFF,1,2,3=ON,1=X,2=Y,3=Z;|| 1=PRESSURE,2=BODY FORCE"<<endl;
out<<"-1 0.0 0.0 0.0		:ACT AT TIMESTEP n, PRESSURE N, PRESSURE P, BODYFORCE(1)"<<endl;
out<<"-1 0.0 0.0 0.0		:ACT AT TIMESTEP n, PRESSURE N, PRESSURE P, BODYFORCE(2)"<<endl;
out<<"-1 0.0 0.0 0.0		:ACT AT TIMESTEP n, PRESSURE N, PRESSURE P, BODYFORCE(3)"<<endl;
out<<"-1 0.0 0.0 0.0		:ACT AT TIMESTEP n, PRESSURE N, PRESSURE P, BODYFORCE(4)"<<endl;
out<<"-1 0.0 0.0 0.0		:ACT AT TIMESTEP n, PRESSURE N, PRESSURE P, BODYFORCE(5)"<<endl;
out<<"0			:BODY FORCE APPLIED ON 1:Phase1,0:BOTH,-1:Phase2"<<endl;
out<<"=====================CAPILLARY PRESSURE MEASURMENTS================"<<endl;
out<<"0 1 3000                     :0=OFF,1,2,3=ON,1=X,2=Y,3=Z;|| 1=PRESSURE,2=BODY FORCE; time steps of 1 changement"<<endl;
out<<"1.0 0.97 9 6           :PRESSURE N,P; Changement times, terminal: no. output intervals "<<endl;                  
out<<"1.0e-4                          :Error of Saturation stable condition"<<endl;

out<<i<<endl;

}


//======================SCRIPT WRITTING========================================
for (int i=1;i<=total_number;i++)
{
        
        name.str("");
	name<<poreFileName2<<i<<".sc";
	ofstream out;
	out.open(name.str().c_str());
	
out<<"#! /bin/csh"<<endl;
out<<"#$ -j y -o $JOB_NAME.$JOB_ID.log"<<endl;
out<<"#$ -cwd"<<endl;
out<<"#$ -pe smp 4"<<endl;
out<<"#$ -l arch=lx24-amd64"<<endl;




out<<"source /usr/local/bin/ifort64vars.csh"<<endl;

out<<"source /usr/local/bin/ompi64vars.csh"<<endl;
out<<""<<endl;

//out<<"cp /work/jy810/Sandpack/LV60/*_"<<i*20000<<"*.bin_input ."<<endl;

out<<endl;

out<<"mpirun -np $NSLOTS paratest "<<poreFileName<<i<<".inputdat";
out<<">report"<<i<<".txt"<<endl;


}


}


